(******************************************************************************)
(**)                        ОТДЕЛ КомплВект;
(******************************************************************************
 * НАЗНАЧЕНИЕ: определение векторов (рядов) комплексных чисел
 *             и разложение их в комплексный ряд Фурье
 *
 * Источники, ссылки, библиография
   Debord J., Библиотека математических подпрограмм
   tpmath.zip

   Goffe B., Программа SimAnn.FOR (Simulated Annealing in Fortran)
   http://www.netlib.org/opt/simann.f

   EISPACK: Библиотека Фортран функций для вычисления собственных значений и векторов
   http://www.netlib.org/eispack

   Marsaglia G., Тесты для генераторов случайных чисел
   http://stat.fsu.edu/~geo/diehard.html

   Moshier S., Библиотека математических подпрограмм
   http://www.moshier.net

   Press W.H., Teukolsky S.A., Vetterling W.T., Flannery B.P.
   Численные рецепты на Паскале
   http://garbo.uwasa.fi/pc/turbopas.html

   Пакет численных методов для Турбо Паскаля
   Borland International, 1986
 *
 ******************************************************************************)
ИСПОЛЬЗУЕТ
  Матем,
  Компл;

ВИД
  Вид-    = РЯД ИЗ Компл.Вид; (* вектор комплексных чисел *)
  Доступ- = ДОСТУП К Вид;
  Вещ     = Матем.Вещ;
  
                              (* заготовки ЗАДАЧ *)
ЗАДАЧА^ БПФ-(вход-,выход+:Вид);
ЗАДАЧА^ ОБПФ-(вход-,выход+:Вид);
ЗАДАЧА^ ПФ-(гармоника:ЦЕЛ; вход-:Вид; выход+:Компл.Вид);

(******************************************************************************)
ЗАДАЧА степень2-(x:ЦЕЛ):ЦЕЛ;
(* Цель: найти такую <степень>, что 2<степень> = < x >
 ******************************************************************************
 * Ответ: <степень> или -1, если её не существует *)
ПЕР
  степень:ЦЕЛ;
УКАЗ
  степень:=0;
  ПОКА ЧЕТ(x) ВЫП
    x:=x ДЕЛИТЬ 2;
    УВЕЛИЧИТЬ(степень)
  КОН;
  ЕСЛИ x # 1 ТО 
    степень:=-1 
  КОН;
  ВОЗВРАТ степень
КОН степень2;

(******************************************************************************)
ЗАДАЧА БитыНаоборот-(биты,всегоБит:ЦЕЛ):ЦЕЛ;
ПЕР
  ц,ответ:ЦЕЛ;
УКАЗ
  ответ:=0;
  ОТ ц:=0 ДО всегоБит-1 ВЫП
    ответ:=(ответ * 2)+(биты ОСТАТОК 2);
    биты:=биты ДЕЛИТЬ 2
  КОН;
  ВОЗВРАТ ответ
КОН БитыНаоборот;

(******************************************************************************)
ЗАДАЧА Преобразовать-(угол:Вещ; вход-,выход+:Вид);
ПЕР
  размерУч,конецУч,размер,всегоБит:ЦЕЛ;
  i,j:ЦЕЛ;
  альфа,бета,прирУгла:Вещ;
  t,арг:Компл.Вид;
УКАЗ
  размер:=РАЗМЕР(вход);
  ПРОВЕРИТЬ(размер > 1);
  всегоБит:=степень2(размер);
  ПРОВЕРИТЬ(всегоБит > 0);
  ОТ i:=0 ДО размер-1 ВЫП
    j:=БитыНаоборот(i,всегоБит);
    (* прореживание входа *)
    выход[j]:=вход[i]
  КОН;
  конецУч:=1;
  размерУч:=2;
  ПОКА размерУч <= размер ВЫП
    прирУгла:=угол/размерУч;
    альфа:=2*Матем.кв(Матем.sin(0.5*прирУгла));
    бета:=Матем.sin(прирУгла);
    i:=0;
    ПОКА i < размер ВЫП 
      Компл.алгВид(1,0,арг); (* (cos(0);sin(0)) *)
      ОТ j:=i ДО i+конецУч-1 ВЫП
        Компл.умн(арг,выход[j+конецУч],t);
        Компл.выч(выход[j],t,выход[j+конецУч]);
        Компл.доб(выход[j],t,выход[j]);
        Компл.алгВид(арг.x-(альфа*арг.x+бета*арг.y),
                     арг.y-(альфа*арг.y-бета*арг.x),арг)
      КОН;
      УВЕЛИЧИТЬ(i,размерУч)
    КОН;
    конецУч:=размерУч;
    размерУч:=размерУч*2
  КОН
КОН Преобразовать;

(******************************************************************************)
ЗАДАЧА БПФ-(вход-,выход+:Вид);
(* Цель: прямое Быстрое Преобразование Фурье (БПФ)
 ******************************************************************************
 * До:    <вход>  - временная последовательность
 * После: <выход> - частотный спектр *)
УКАЗ
  Преобразовать(Матем.ПИ2,вход,выход)
КОН БПФ;

(******************************************************************************)
ЗАДАЧА ОБПФ-(вход-,выход+:Вид);
(* Цель: обратное Быстрое Преобразование Фурье (оБПФ)
 ******************************************************************************
 * До:    <вход>  - частотный спектр
 * После: <выход> - временная последовательность *)
ПЕР
  размер,i:ЦЕЛ;
УКАЗ
  Преобразовать(-Матем.ПИ2,вход,выход);
  (* нормировка полученного временного ряда *)
  размер:=РАЗМЕР(вход);
  ОТ i:=0 ДО размер-1 ВЫП
    Компл.делВ(выход[i],размер,выход[i])
  КОН
КОН ОБПФ;

(******************************************************************************)
ЗАДАЧА ПФ-(гармоника:ЦЕЛ; вход-:Вид; выход+:Компл.Вид);
(* Цель: Преобразование Фурье для одной заданной гармоники
 ******************************************************************************
 * До:    <гармоника> - номер гармоники
 *        <вход>      - временная последовательность
 * После: <выход> - коэффициенты разложения на заданной гармонике *)
ПЕР
  тета,бета:Вещ;
  размер,i:ЦЕЛ;
  к1,к2,к:Компл.Вид;
УКАЗ
  размер:=РАЗМЕР(вход);
  выход:=Компл.НОЛЬ;
  тета:=Матем.ПИ2*гармоника/размер;
  Компл.алгВид(Матем.cos(-2*тета),Матем.sin(-2*тета),к1);
  Компл.алгВид(Матем.cos(-тета),  Матем.sin(-тета)  ,к2);
  бета:=2*к2.x;
  ОТ i:=0 ДО размер-1 ВЫП
    Компл.умнВ(к2,бета,к);
    Компл.выч(к,к1,к);
    к1:=к2; к2:=к;
    Компл.умн(вход[i],к,к);
    Компл.доб(выход,к,выход)
  КОН
КОН ПФ;

КОН КомплВект.